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In this talk we review the GOLEM approach to one-loop calculations and present an automated 
implementation of this technique. This method is based on Feynman diagrams and an advanced 
reduction of one-loop tensor integrals which avoids numerical instabilities. We have extended 
our one-loop integral library golem 9 5 with an automated one-loop matrix element generator 
to compute the virtual corrections of the process qq — > bbbb. The implementation of the virtual 
matrix element has been interfaced with tree-level Monte Carlo programs to provide the full result 
for the above process. 
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1. Overview 

The successful start of the LHC and its first data taking mark the dawn of a new era in particle 
physics. Very soon we will be able to explore the energy scale of electro-weak symmetry breaking 
which enables us to either confirm the Standard Model as a low energy theory of particle physics 
or to discover new particles guiding us to extensions of the Standard Model. No matter what 
the outcome of the experiment will be, a successful interpretation of the data will require a large 
number of predictions calculated at least to Next-To-Leading Order (NLO) in QCD. Some of these 
reactions involve up to four particles in the final state [jlj |2|, ^, |p. 

The computation of NLO corrections to processes involving many final state particles are 
cumbersome and time consuming. The time required for a single calculation — without the use of 
automated programs — often coincides with the duration of a Ph.D. Considering the fast experi- 
mental progress it is very likely that already within the next two or three years the physics results 
of the LHC will raise the demand for precision calculations within and beyond the Standard Model 
for many multileg processes. Only the automatisation of NLO calculations will allow us to keep 
pace with the requirements set by the experiments. 

We argue that the automatisation of the computation of NLO cross-sections will also improve 
the possibility of comparing results from different implementations, especially if the tools are made 
public and common conventions are used for input/output and interfacing to external programs are 
in use. A first step in this direction has been made through the Binoth Les Houches Accord [0]. 
This accord exploits the modular structure of NLO calculations and proposes the reflection of this 
structure in the implementation of such calculations as computer programs. 

Any QCD cross-section at NLO can be written in the form 

The first term on the right-hand side describes the Born-level cross-section calculated from the 
squared tree-level amplitude of the 2 — > N process. The virtual corrections o{^ N stem from the 
interference term between tree-level and one-loop diagrams. The third term describes the real 
radiation of an extra, unobserved parton at tree-level. 

The last two terms lead to singularities which can be regularised by introducing a non-integer 
dimension n = 4 — 2e, yielding poles in l/e which cancel only after both terms have been added 
up. For practical applications it is therefore convenient to introduce subtraction terms, such that 
both o£+ N an d Oz^N+i are fi mte an d can be integrated over phase space independently. 

A full implementation of an NLO calculation can therefore be modularised into a Monte- 
Carlo integrator for phase space integration, a tree-level matrix element generator for of_^ N and 
of^ N+l , a one-loop matrix element generator for o^ N and infrared subtraction terms. Typically 
the first two of these components are implemented in the same program. An overview of existing 
techniques and recent contributions to all of the modules can be found, for example, in [Q] and [j|]. 

Methods for computing o^n can ^ e classified in two groups. On the one hand there are 
unitarity based methods determining the coefficients of scalar one-loop integrals by exploiting 
analyticity of the amplitude; on the other hand there are Feynman diagrammatic techniques starting 
from tensor integrals, which are then reduced to simpler integrals which can be evaluated in a 
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numerically stable way. Both techniques have led to automated implementations in recent years 
([§, 0, g g [1(1 Ol and [@ HH [If]] respectively). 

In Sections ^| and [| we present the automated implementation of a method using Feynman di- 
agrams for the calculation of Oj-jjy, tne GOLEM method, based on the reduction scheme proposed 
in [17]. This implementation has been used to compute the NLO virtual corrections to qq — > bbbb, 
which is one of the two partonic initial states contributing to the pp — > bbbb process. This pro- 
cess constitutes an important background for Higgs searches in models beyond the Standard Model 
where Higgs bosons decay predominantly into ^-quarks as discussed in [ 1 8 , II 



2. The Generic Integral Form Factor Library golem95 

The Feynman diagrammatic approach to computation of one-loop corrections to processes 
with N external particles requires the evaluation of tensor integrals of rank r which have the general 
form 

u ft qa\ qa r 

i% n l 2 {q\ -m\ + id)---{q 2 H -m] i + id) ' 

where qt = k + r; and S denotes the matrix 5,; = (r,- — rj) 2 — mf — m 2 . 

It is well-known [ |20| ] that such a tensor integral can be expressed in terms of basis scalar inte- 
grals, but the reduction procedure introduces inverse Gram determinants (detG) in the coefficients 
of the expansions which can lead to numerical instabilities in certain regions of the phase space. 



7 N (a u ...,a r ,S)- I -—-—-—————j ■ — , (2.1) 



Therefore, in Q21| , |17| ] we have proposed a reduction scheme which allows to write any N- 
point amplitude as a linear combination of basis integrals (I2, 1", 7" +2 , 7^ +2 , 1'l +A ) with and without 
Feynman parameters in the numerator, avoiding the introduction of inverse Gram determinants. 

The evaluation of the basis functions can be performed by reducing them further to scalar 
integrals using recursion formulae. This further reduction, however, introduces Gram determinants 
in the coefficients which could lead to numerical instabilities in certain regions of the phase space. 

Potentially dangerous regions are identified by the criterion | detG| < A| det»S| for a fixed cut- 
off A ~ 10~ 5 . In this regions the basis integrals are evaluated numerically without applying any 
further reduction. 

For all integrals without internal masses we have worked out one-dimensional integral repre- 
sentations which can be evaluated by numerical integration. 

This algorithm has been implemented in form of a Fortran 90 library, golem 9 5, for massless 
internal propagators (m, = 0) [jl^]. This version of the code has been made available for download 1 . 

We have recently extended the library golem9 5 to the case where internal masses are present. 
All infrared divergent integrals have been implemented explicitly. For evaluating the finite boxes 



and triangles the user needs to link the LoopTools library [ |12| , [22| , [23| ]. This "massive" version of 
the golem 9 5 library is currently in the testing phase and will be available shortly. 

For integrals with internal masses, the option to evaluate the tensor integrals numerically prior 
to reduction in regions where the Gram determinant tends to zero, is not yet supported. How- 
ever, one-dimensional integral representations valid for all possible kinematic configurations are 
under construction. 



http: //lappweb. in2p3 . f r /lapth/GOLEM/golem95 . html 
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3. Implementation of a One-Loop Matrix Element Generator 

Building upon golem9 5 as a loop-integral library, our next step was the construction of a 
matrix-element generator at the one-loop level. The computation is carried out projecting the am- 
plitude onto helicity and colour structures. The virtual corrections can therefore be expressed as 

do£*=— £ ^(^,^^h...,^)^c y |c,)^(^,p^pfS...,^)+^c. (3.1) 
n " nb {k},j,k 

where p i ' denotes the pair of momentum /?, and helicity label A, of the i— th particle. The matrix 
(cj\ck) consists of the contractions of all colour basis tensors for a given process evaluating to 
rational functions in the number of colours Nc and the normalisation constant Tr of the generators. 
The constants n a and nj, represent the averaging over spin and colour. 

The one-loop amplitude Y*k^k\ck) consists of a sum of Feynman diagrams, which we generate 
using QGraf []24"|]. We use QGraf together with LTjjX and Axodraw [[25|] also for drawing the 
diagrams; the layout of the diagrams is determined using the algorithm of 



The expressions of the diagrams are then processed using Form [ |27| ] and the Form library 
spinney which we have developed for dealing with helicity spinors and ^-dimensional Dirac and 
Lorentz algebra efficiently. Majorana spinors can also be dealt with thanks to the the implementa- 
tion the flipping rules for spin lines as described in [|28|]. 

At the moment the GOLEM program can import CompHep [ ]29| ] model files to perform Beyond 
the Standard Model computations. An interface to FeynRules [^(J is under construction. 

After the Form program has decomposed the diagram expression into colour structures and 
the tensor integrals are represented in terms of integral form factors as defined in golem95 [|TJ], 
the resulting expressions are optimised and translated into Fortran 90 functions using the code 



generator haggles pip. At this step the number of multiplications is minimised applying a 
Horner scheme and common subexpression elimination. 

The generated Fortran 90 program is linked with go lem9 5 for the numerical evaluation of the 
tensor integral form factors. A future version of the program will support the Binoth Les Houches 
Accord []|] to facilitate the interfacing to Monte-Carlo generators. 

4. NLO Results for qq -»■ bbbb [{32p 

The setup described in the previous section has been used to compute the virtual corrections 
of the QCD corrections to qq — > bbbb in the limit nit = and m t — > °°. We have compared the 
results with an independent implementation using FeynArts and FormCalc [|l^] to generate and 
simplify the diagrams, where the tensor integrals are algebraically reduced to scalar integrals using 



the procedure described in []17[]. 

In order to compute qq — > bbbb at NLO accuracy the Born level cross-section, the real emis- 
sion contribution and the infrared subtraction terms also need to be evaluated. Since we are only 
interested in a process with 4 tagged Z?-jets the relevant process for the real emission contribution 
is qq — > bbbbg. 

We have used MadGraph/MadEvent [33|, 34] and MadDipole p"5p to evaluate the tree- 
like contributions and to perform the phase space integration. As an alternative setup, based on 
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an extended version of the Whizard [3f| 37] Monte Carlo generator with an implementation of 
infrared subtraction terms has been used to obtain an independent cross-check. For the infrared 
subtraction we have used Catani-Seymour dipoles [^] in both implementations including a phase 
space slicing parameter following 



To define a bbbb event, we first we apply a kj jet algorithm [ 10] to decide if the gluon should 
be merged with a &-quark. If the gluon is merged we use the effective momentum py = py + p g 
as the momentum of the ft-quark in the cuts and in the observables. Then we apply a pj cut of 
Pr(bj) > 30GeV and a rapidity cut of \rj \ < 2.5 to all ft-quarks and a separation cut of AT? > 0.8 to 
all pairs of ft-quarks. 



We sum over q £ {u,d,c,s} and use the CTEQ6M parton distribution functions plQ with two- 
loop running of CC S for both Leading and Next-to-Leading Order computations. The centre of mass 
energy is set to y/s = \4GeV. In our results we use a fixed factorisation scale of jip = lOOGeV; 



the renormalisation scale we set to /lo = y Ly Pj{bj). 

In Figure [j] we show the invariant mass distribution of the system of the two Z?-pairs with 
highest pj. The error bands have been obtained by varying the renormalisation scale }Xr = x/Xq 
between 1/4 < x < 2. The dashed line marks the leading order distribution for x = 1/2, which 
turns out to be very similar to the NLO prediction for this value. The reduction of the uncertainty 
band due to scale variations clearly shows the importance of the NLO corrections for the precision 
of this calculation. 



> 



E 




350 400 
m(b i ,b 2 ) [GeV] 



Figure 1: Invariant mass {m^) disnibution of the two leading b-quarks (see text). The error bands are 
obtained by varying the renormalisation scale /Xr = x/Xo between 1/4 < x < 2, where jtio = Pr^Pj)- The 
dashed line shows the value of the leading order prediction for x =1/2. 
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5. Conclusion 

We have presented results obtained using the GOLEM method and recent progress in the im- 
plementation of an automated one-loop matrix element generator. A first important step towards 
this goal was the development of a one-loop integral library, golem9 5, which is currently being 
extended to the case of massive propagators. As a second step we have added a completely auto- 
mated framework which generates efficient Fortran 90 code for the numerically stable evaluation of 
one-loop matrix elements from a set of Feynman rules. This framework includes the development 
of new tools such as an optimising code generator (haggies), and a Form library (spinney) for 
the ueatment of helicity spinors and 7i-dimensional Dirac and Lorentz algebra. Finally, we have 
presented the complete NLO result for the process qq — > bbbb, which is a subprocess of — >• bbbb, 
an important background to Higgs searches beyond the Standard Model. 

In the near future we plan to implement the interface to Monte-Carlo tools described in [||] and 
to make all parts of the program publicly available. 
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